Sun & Breheny (2019) webgazer replication

This study is a replication of Sun & Breheny (2019) with the web-based eye-tracking paradigm using the javascript webgazer library.

Here’s the experiment we’re analyzing.

library(tidyverse)
library(lme4)
library(ggplot2)
library("readxl")

# load helper scripts
source("helpers.R")

# load target word onsets
audio_info = read_excel("../data/audiotimes.xlsx")

# set color-blind-friendly color palette
cbPalette <- c("#56B4E9", "#D55E00", "#009E73","#999999", "#E69F00","#009E73","#56B4E9", "#D55E00", "#009E73","#999999", "#E69F00","#009E73","#56B4E9", "#D55E00", "#009E73","#999999", "#E69F00","#009E73","#56B4E9", "#D55E00", "#009E73","#999999", "#E69F00","#009E73")

# load data 
d = read.csv('../data/sunbreheny_webgazer_pilot2-merged.csv') %>%
  select(-location1,-location2,-location3,-location4,-location5,-location6,-location7,-location8,-location9,-location10,-target_figure,-target_object,-target_num_object,-error) %>% 
  filter(trial_type=="Exp") %>%
  mutate(target_location = ifelse((target1==7 | target2==7),"topleft",ifelse((target1==8 | target2==8),"bottomleft",ifelse((target1==9 | target2==9),"topright",ifelse((target1==10 | target2==10),"bottomright",NA))))) %>%
  mutate(competitor_location = ifelse((competitor1==7 | competitor2==7),"topleft",ifelse((competitor1==8 | competitor2==8),"bottomleft",ifelse((competitor1==9 | competitor2==9),"topright",ifelse((competitor1==10 | competitor2==10),"bottomright", NA))))) %>%
  left_join(audio_info,by=c("audio" = "Prime")) %>% #for testing: TODO
  filter(as.character(proliferate.condition)!= "list1") %>%
  filter(as.character(proliferate.condition)!= "list2") %>%
  filter(as.character(proliferate.condition)!= "list3")

# formatting
d$webgazer_time = gsub("\\[","",d$webgazer_time)
d$webgazer_time= gsub("\\]","",d$webgazer_time)
d$x = gsub("\\[","",d$x)
d$x = gsub("\\]","",d$x)
d$y = gsub("\\[","",d$y)
d$y = gsub("\\]","",d$y)

Exclude trials with wrong selections and plot subject information.

d = d %>% 
  mutate(selection_correct = ifelse(as.character(response) == as.character(target1),1, ifelse(as.character(response) == as.character(target2),1,0)),language = tolower(subject_information.language))

table(d$selection_correct)

  0   1 
  1 107 

d = d %>% 
  filter(selection_correct==1)

participants = d %>% 
  select(workerid,subject_information.accuracy,subject_information.age,subject_information.camblock,subject_information.comments, subject_information.eyesight, subject_information.eyesight_task, subject_information.gender, subject_information.headphones, language, subject_information.previous_accuracy_attempts, subject_information.time_in_minutes) %>% 
  unique()

# number of unique participants:
nrow(participants) 

[1] 3

# distribution of participant gender
table(participants$subject_information.gender)

 Male Woman 
    2     1 

# distribution of completion times
summary(participants$subject_information.time_in_minutes)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  12.11   12.82   13.52   13.32   13.92   14.31 

ggplot(participants, aes(x=subject_information.time_in_minutes)) +
  geom_histogram()

# distribution of calibration accuracy
summary(participants$subject_information.accuracy)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  55.00   63.50   72.00   66.67   72.50   73.00 

ggplot(participants, aes(x=subject_information.accuracy)) +
  geom_histogram()

# distribution of ages
summary(participants$subject_information.age)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  24.00   25.00   26.00   40.33   48.50   71.00 

ggplot(participants, aes(x=subject_information.age)) +
  geom_histogram()

# eyesight information
ggplot(participants, aes(x=subject_information.eyesight,fill=subject_information.eyesight_task)) +
  geom_histogram(stat="count")

# distribution of camera views -- some people continued to see camera view even during experiment (bug we haven't figured out yet how to fix)
ggplot(participants, aes(x=subject_information.camblock,fill=subject_information.camblock)) +
  geom_histogram(stat="count")

# distribution of headphone use
ggplot(participants, aes(x=subject_information.headphones)) +
  geom_histogram(stat="count")

# distribution of native languages
ggplot(participants, aes(x=language)) +
  geom_histogram(stat="count") +
  theme(axis.text.x=element_text(angle=45, hjust=1,vjust=1))

Separate time and x/y coordinates into separate data points, one per row, and define regions of interest (ROIs).

Bin samples and align time to the onset of words.

# relevant time columns------------
# trial_start: unix time at the beginning of the trial (when images are displayed)
# audio_play_unix: unix time at the beginning of the audio (after 1sec display preview)
# webgazer_time: elapsed time since webgazer started
# unix_tlist: unix time of each webgazer_time sample
# unix_rt: unix time of selection
# response_time: response time in ms (trial_start - unix_rt)

first_samples = dd %>% 
  group_by(workerid,slide_number) %>% 
  summarize(first_sample = min(webgazer_time))

# total number of samples:
nrow(first_samples)

# set size of time bins to collapse over
binsize = 60

# TODO: FIX EVERYTHING AFTER THIS --- 
# bin samples and align to target word onset
dd_binned = dd %>% 
  left_join(first_samples,by=c("workerid","slide_number")) %>% 
  mutate(relative_time = webgazer_time-first_sample) %>% 
  mutate(relative_time_bin=floor(relative_time/binsize)) %>% 
  mutate(binned_time_relative = relative_time_bin*binsize)

Plot proportions of looks to all regions (target, competitor, center, other)

# aggregate looks (compute proportions of looks to each region in each time bin and condition; add 95% bootstrapped confidence intervals)
agr = dd_binned %>% 
  group_by(condition,binned_time_relative) %>% 
  summarize(target_prop=mean(target_look),target_CILow=ci.low(target_look),target_CIHigh=ci.high(target_look),competitor_prop=mean(competitor_look),competitor_CILow=ci.low(competitor_look),competitor_CIHigh=ci.high(competitor_look),center_prop=mean(center_look),center_CILow=ci.low(center_look),center_CIHigh=ci.high(center_look),other_prop=mean(other_look),other_CILow=ci.low(other_look),other_CIHigh=ci.high(other_look)) %>% 
  ungroup() %>% 
  mutate(target_ymin=target_prop-target_CILow,target_ymax=target_prop+target_CIHigh,competitor_ymin=competitor_prop-competitor_CILow,competitor_ymax=competitor_prop+competitor_CIHigh,center_ymin=center_prop-center_CILow,center_ymax=center_prop+center_CIHigh,other_ymin=other_prop-other_CILow,other_ymax=other_prop+other_CIHigh)

# prepare data for plotting
long_props = agr %>% 
  select(condition,binned_time_relative,target_prop,competitor_prop,center_prop,other_prop) %>% 
  pivot_longer(cols = target_prop:other_prop,names_to=c("region"),values_to=c("proportion")) %>% 
  separate(region,c("region",NA))

long_ymin = agr %>% 
  select(condition,binned_time_relative,target_ymin,competitor_ymin,center_ymin,other_ymin) %>% 
  pivot_longer(cols = target_ymin:other_ymin,names_to=c("region"),values_to=c("ymin")) %>% 
  separate(region,c("region",NA))

long_ymax = agr %>% 
  select(condition,binned_time_relative,target_ymax,competitor_ymax,center_ymax,other_ymax) %>% 
  pivot_longer(cols = target_ymax:other_ymax,names_to=c("region"),values_to=c("ymax")) %>% 
  separate(region,c("region",NA))

toplot = long_props %>% 
  left_join(long_ymin,by=c("condition","binned_time_relative","region")) %>% 
  left_join(long_ymax,by=c("condition","binned_time_relative","region")) %>% 
  mutate(region = fct_relevel(region,"target","competitor","center"))